Plot the Dow Jones daily by using R built-in functions:
##
## Attaching package: 'psych'
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
df = 3
tmp = ggplot(DJ_d, aes(sample=r_Dow_Jones/sd(r_Dow_Jones)*sqrt(df/(df-2))))
tmp + geom_qq(distribution=stats::qt, dparams=list(df=df)) + geom_abline(aes(intercept=0,slope=1),color="red")
#qqplot(rt(length(DJ_d$r_Dow_Jones),df=5),DJ_d$r_Dow_Jones)
vq = pnorm(vdata)
hist(vq)
ggplot(data=tibble(vq),aes(vq))+geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
vn = NULL; for(val in grids) vn = c(vn,sum(vq <= val))
vn = c(vn[1],diff(vn))
test = sum((vn-length(vdata)/ik)**2/(length(vdata)/ik))
cat("test =",test," df =",ik-3," p-value =",1-pchisq(test,df=ik-3))
## test = 3405.87 df = 17 p-value = 0
df = 3
ndata = vdata*sqrt(df/(df-2))
vq = pt(ndata,df=5); hist(vq)
vn = NULL; for(val in grids) vn = c(vn,sum(vq <= val))
vn = c(vn[1],diff(vn))
test = sum((vn-length(vdata)/ik)**2/(length(vdata)/ik))
cat("test =",test," df =",ik-3," p-value =",1-pchisq(test,df=ik-3))
## test = 414.7555 df = 17 p-value = 0
alpha = 0.1514736
sigma = 4.0013995
ndata = vdata * sqrt((1-alpha) + alpha*sigma**2)
vq = (1-alpha)*pnorm(ndata) + alpha*pnorm(ndata,sd=sigma); hist(vq)
vn = NULL; for(val in grids) vn = c(vn,sum(vq <= val))
vn = c(vn[1],diff(vn))
test = sum((vn-length(vdata)/ik)**2/(length(vdata)/ik))
cat("test =",test," df =",ik-3," p-value =",1-pchisq(test,df=ik-3))
## test = 212.4475 df = 17 p-value = 0
Write a function for minimization to choose the optimal pair for \(\alpha\) and \(\sigma\).
func <- function(vx){
alpha = vx[1]
sigma = vx[2]
ndata = vdata*sqrt((1-alpha) + alpha*sigma**2)
vq = (1-alpha)*pnorm(ndata) + alpha*pnorm(ndata,sd=sigma)
vn = NULL; for(val in grids) vn = c(vn,sum(vq <= val))
vn = c(vn[1],diff(vn))
return(sum((vn-length(vdata)/ik)**2/(length(vdata)/ik)))
}
func(c(0.15,4))
## [1] 212.7936
optim(par=c(0.15,4),fn=func,method="BFGS")
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## Warning in sqrt((1 - alpha) + alpha * sigma^2): NaNs produced
## $par
## [1] 0.1515383 4.0008180
##
## $value
## [1] 212.4475
##
## $counts
## function gradient
## 54 2
##
## $convergence
## [1] 0
##
## $message
## NULL
optim(par=c(0.15,4),fn=func,method="L-BFGS-B",lower=c(0,0),upper=c(1,10))
## $par
## [1] 0.1752619 3.7378776
##
## $value
## [1] 190.9963
##
## $counts
## function gradient
## 148 148
##
## $convergence
## [1] 52
##
## $message
## [1] "ERROR: ABNORMAL_TERMINATION_IN_LNSRCH"
atmp = seq(from=0.05, to=.95, length.out=40)
stmp = seq(from=1, to=10, length.out=50)
tmp = expand.grid(atmp,stmp)
res = apply(tmp,1,func)
res = t(matrix(res, nrow=length(atmp)))
library(plotly)
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
ret = plot_ly(x=atmp, y=stmp, z=res) %>% add_surface() %>%
layout(scene=list(xaxis=list(title='alpha'), yaxis=list(title='sigma'), zaxis=list(title='function')))
ret
lret = apply(log(data),2,diff)
summary(lret)
## DAXINDX FRCAC40 FTSE100
## Min. :-0.1370990 Min. :-0.1014 Min. :-0.1302860
## 1st Qu.:-0.0053976 1st Qu.:-0.0060 1st Qu.:-0.0045157
## Median : 0.0002581 Median : 0.0000 Median : 0.0002738
## Mean : 0.0003493 Mean : 0.0003 Mean : 0.0004102
## 3rd Qu.: 0.0070188 3rd Qu.: 0.0069 3rd Qu.: 0.0058174
## Max. : 0.0728754 Max. : 0.0823 Max. : 0.0759697
## NA's :393
## HNGKNGI NIKKEI SNGALLS
## Min. :-0.4054220 Min. :-1.614e-01 Min. :-0.0940305
## 1st Qu.:-0.0053347 1st Qu.:-6.014e-03 1st Qu.:-0.0041041
## Median : 0.0002220 Median : 0.000e+00 Median : 0.0000000
## Mean : 0.0005713 Mean : 4.991e-05 Mean : 0.0001921
## 3rd Qu.: 0.0079703 3rd Qu.: 6.380e-03 3rd Qu.: 0.0047405
## Max. : 0.1724710 Max. : 1.243e-01 Max. : 0.1431300
##
## SPCOMP
## Min. :-0.2283303
## 1st Qu.:-0.0033884
## Median : 0.0003843
## Mean : 0.0004885
## 3rd Qu.: 0.0049146
## Max. : 0.0870888
##
matplot(lret,type='l',ylab='log returns',xlab='time horizon')
sum(is.na(lret[,'FRCAC40']))
## [1] 393
LB <- function(vx,lag,ip){
tmp = acf(vx,lag.max=lag,plot=F)$acf
tmp = tmp[2:(lag+1)]**2
test = sum(tmp/(length(vx)-1:lag))*length(vx)*(length(vx)+2)
return(list(test=test, pval=1-pchisq(test,df=lag-ip)))
}
tmp = 6
LB(vx=lret[!is.na(lret[,tmp]),tmp],lag=100,ip=0)
## $test
## [1] 202.3019
##
## $pval
## [1] 6.465107e-09